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Abstract. This is part two in a series of papers in which we investigate an approach based on 
£S| Lagrangian perturbation theory (LPT) to study the non-linear evolution of the large-scale 

structure distribution in the universe. Firstly, we compute the matter bispectrum in real 

space using LPT up one-loop order, for both Gaussian and non-Gaussian initial conditions. 
{N) In the initial position limit, we find that the one-loop bispectrum computed in this manner 

is identical to its counterpart obtained from standard Eulerian perturbation theory (SPT). 

Furthermore, the LPT formalism allows for a simple reorganisation of the perturbative series 

corresponding to the resummation of an infinite series of perturbations in SPT. Applying 
j—{ this method, we find a resummed one-loop bispectrum that compares favourably with results 

from A-body simulations. We generalise the resummation method also to the computation 
•'"I of the redshift-space bispectrum up to one loop. 
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1 Introduction 

Current measurements of the cosmic microwave background (CMB) anisotropies and large- 
scale structure (LSS) distribution strongly support the validity of the so-called ACDM 
model [1]. In this model, quantum fluctuations on an inflaton field set the seeds of pri- 
mordial curvature perturbations in the very early universe. These perturbations manifest 
themselves as inhomogeneities in the matter/energy fields, and are observed as temperature 
and polarisation anisotropies in the CMB. Subsequent evolution via gravitational instability 



significantly enhances these initial perturbations, ultimately leading to the formation of the 
cosmic structures we see today. 

Although the basic ACDM paradigm has so far been remarkably successful at explaining 
a host of astrophysical observations, it is nonetheless crucial to devise more tests to constrain 
the model's parameters, to find its boundaries of validity, and to obtain more insight into the 
details of its sub-scenarios. One interesting issue is whether or not the primordial curvature 
perturbations conform to a perfect Gaussian distribution. (Almost) perfect Gaussianity is 
a hallmark feature of a large class of simple inflation models, namely, the canonical single- 
field slow-roll models [2]. Nonetheless, many other observationally consistent inflation models 
are capable of producing primordial perturbations that deviate significantly from Gaussianity 
(i.e., primordial non-Gaussianity; PNG) [3]. Therefore, any future measurement or constraint 
of PNG would strongly limit the inflationary model space. 

In general, the term "Gaussianity" for a random variable 5(x) can be defined in terms 
of the distribution of complex phases <f> from the Fourier decomposition of an ensemble of 
random realisations. A uniform distribution in [0, 2tt] indicates a completely randomised, 
and hence Gaussian, distribution; any deviation from uniformity is to be identified with 
"non-Gaussianity" [4]. The much- used two-point statistics, the power spectrum — defined as 
P(k) = \5(k)\ 2 , where k = \k\ is the Fourier wavenumber — is not well-suited to study non- 
Gaussianities, because it is independent of the phases of the perturbation variable 5(k). The 
lowest order statistics that is directly sensitive to the phases is the three-point function, or the 
bispectrum £?(fa, fa, fa), defined via (S(ki)6(k2)S(k3)) c = (2ir) 3 5 D ' (ki+k2+ks)B(ki, fa, fa). 
Gaussian fluctuations necessarily lead to B(kx, fa, fa) = 0, since the sum of the phases <£ = 
<p(k±) + </>(fa) + </>(fa) must also be uniformly distributed in [0, 27r], so that the ensemble 
average (e ) evaluates to zero. Contrastingly, a non-uniform, and hence non-Gaussian, 
phase distribution always yields a non-zero -B(fa, fa, fa). 

Although a non-zero bispectrum necessarily indicates the presence of non-Gaussianity, 
it is important to realise that non-Gaussianities arise generically as a result of non-linear 
coupling between different Fourier modes. This means that while linear evolution of the 
primordial curvature perturbations preserves their statistical properties, as soon as linear 
theory fails to describe the evolution of inhomogeneities, a sizable bispectrum can be expected 
from non-linear evolution alone, even in the absence of PNG. A particularly relevant case 
is the clustering statistics of the present-day LSS distribution; non-linear evolution of the 
matter density perturbations at low redshifts necessarily generate non-Gaussianities of its 
own. Therefore, if we wish to use the LSS bispectrum as a probe of PNG, it is important 
that we know how to filter out this "late-time" contribution. 

Fully non-linear evolution of the matter density perturbations is in general not amenable 
to analytical treatments. However, if we are merely interested in the mildly non-linear 
regime — where linear theory still dominates, and non-linear evolution contributes a small 
correction to the linear solution — then semi-analytic techniques based on solving a set of 
fluid equations using a perturbative expansion generally return reasonable results (e.g., [5, 
6]). Non- linear corrections to the LSS power spectrum up to two loops [7-9] have been 
computed within the framework of standard perturbation theory (SPT) in, e.g., [10-14], while 
references [15, 16] work within the framework of Lagrangian perturbation theory (LPT). The 
bispectrum, on the other hand, has only been treated using SPT up to one loop [17, 18]. 

In this work, we compute for the first time the tree-level contribution and the one-loop 
correction to the LSS bispectrum using the framework of Lagrangian perturbation theory. 
We use the perturbative solutions up to fourth order from our accompanying paper [19], and 



show that, in the so-called initial position limit, the one-loop bispectrum obtained using the 
LPT formalism is identical to its SPT counterpart. Furthermore, an advantage of the LPT 
formalism is that the perturbative series can be easily reorganised, so that an infinite series of 
perturbations in SPT is effectively resummed [15]. Generalising this resummation method to 
the bispectrum calculation, we find a "resummed" one- loop bispectrum that, when compared 
with a naive expansion, is generally a better approximation to the "exact" bispectrum ex- 
tracted from TV-body simulations in the mildly non-linear regime. Resummation techniques 
have been explored in many matter power spectrum calculations [5, 6, 13-16, 20]. To our 
knowledge, they have as yet not been applied to the computation of the bispectrum. 

The paper is organised as follows. In section 2 we review the formalism of LPT, and write 
down the general expression for the bispectrum. We evaluate this expression in section 3 up 
to one-loop order, using both a simple Taylor expansion and the aforementioned resummation 
technique. The resulting LPT bispectrum and resummed bispectrum are then compared with 
data from iV-body simulations. In section 4, we generalise the resummation technique to the 
computation of the redshift-space bispectrum up to one-loop order. Our conclusions can be 
found in section 5. In general we try to keep the technical details in the main text to the 
minimum necessary for the sake of readability; the reader will be referred to the appendices 
for the details of the computations at the appropriate points. Here we highlight especially 
appendix A, where we report all perturbative kernels in both LPT and SPT up to fourth 
order. 

2 Formalism 

We briefly review in this section the concepts of Lagrangian perturbation theory (LPT), and 
construct the bispectrum using the central object of LPT, the displacement field \l/. Since 
the scales that require non-linear corrections are generally well inside the Hubble horizon, 
a full general relativistic treatment is not necessary, and we work within the Newtonian 
limit of cosmological perturbation theory [21]. For readers wishing to skip to the crux of this 
work, the starting base of our bispectrum calculation is equation (2.17), while equation (2.22) 
shows the reorganised perturbative series that is the starting expression for the resummed 
bispectrum. 

2.1 Lagrangian perturbation theory 

In the Lagrangian framework, the observer follows the trajectories of the individual fluid 
elements [22], where each trajectory is encoded in the time-integrated displacement field \JJ. 
The comoving position a; of a fluid element at conformal time r is then given by its initial 
Lagrangian coordinate q plus the displacement field \l/ evaluated at the same time: 

x(q,T)=q + V(q,T). (2.1) 

The volume of the specific fluid element generally deforms as a result of gravitationally 
induced displacement. The deformation is encoded in d 3 x = J(q, r) d 3 q, where the Jacobian 
determinant J = det[dx/dq] is generally a function of \1/. Mass conservation (for a non- 
relativistic fluid) then leads to the constraint: 

p(x,r)d 3 x = p(T)d 3 q, (2.2) 

where p(x,r) = p(t) [1 + 5(x,t)] is the Eulerian density field with density contrast 5(x,t), 
and p(r) is the mean mass density. Note that in writing down the constraint (2.2), we 



have assumed the initial density contrast S(q, Tmi) at q to be negligibly small [19]. Simple 
algebraic manipulations of equation (2.2) then allow us to relate the density contrast to the 
displacement field: 

In the following, where there is no confusion, we will omit writing out the time dependence 
of the dynamical variables. 

It remains to write down an equation of motion for \l/. In a non-rotating (Eulerian) 
frame, the equations of motion for self-gravitating and irrotational dust are [23-25]: 

^ + V x -{[1 + 5(x,t)]u(x,t)} = 0, (2.4) 



dw „., \ , x 

— + H{t)u(x,t) 



3 -,2 



2 H'6(x,t), (2.5) 



V x xu{x,t) = 0, (2.6) 



where u = dx/dr is the peculiar velocity, d/dr = d/dr+u-'V x the convective derivative, and 
the conformal Hubble parameter evaluates to % = 2/r in a matter-dominated flat universe. 
The last equation (2.6) states the irrotationality of the fluid motion in Eulerian space. We 
call equations (2.4) to (2.6) the Euler-Poisson system. 

To turn the Euler-Poisson system into equations of motion for \I/, we first convert the 
Eulerian derivative Va, to its Lagrangian counterpart using the Jacobian of the coordinate 
transform, while formally requiring it to be invertible. 1 Defining Jij := Sij + VPjj, where the 
subscript ", j" denotes a partial derivative with respect to the Lagrangian coordinate qj, the 
resulting set of equations contains only one dynamical variable, namely, the displacement 
field * [261: 



d 2 J jti (V) _ 12 
d.77 2 r/ 2 



eu m e jpq J Pt i(*)Jq, m (*) f\ ' = -5 [J(¥) - 1] , (2.7) 



A J. ,(qr) 
Ji,n(*) £njkJlj(*) d = , with J > , (2.8) 

where drj = dr/a is the super conformal time, and summation over repeated indices is implied. 
We call this closed set of equations the Lagrange-Newton system. 

Importantly, the Eulerian irrotationality condition (2.6) gives rise to a set of Lagrangian 
constraints on the displacement field, i.e., equation (2.8); in general, \l/ contains both a lon- 
gitudinal and a transverse component, that latter arising from the non-linear transformation 
between the Eulerian and the Lagrangian frame. Lagrangian transverse fields are a necessary 
constituent of the total displacement field, as they ensure irrotationality of the fluid motion 
in the Eulerian frame. This means that both (Lagrangian) longitudinal and (Lagrangian) 
transverse parts affect equation (2.7), and hence the (Eulerian) longitudinal part of the Euler 
equation (2.4) as well. Thus, there is in general no decoupling between the longitudinal and 
the transverse components. 

Equations (2.7) and (2.8) can be solved using a perturbative series \1/ = Yl^Li ^ • The 
linear (n = l) solution in the initial position limit 2 is simply the Zel'dovich approximation [19, 



1 Lagrangian solutions exist even if the Jacobian is not invertible. However, equivalence between the Euler- 
Poisson system and the Lagrange-Newton system can be established only if J 7^ [19]. 

2 The initial position limit corresponds to the case in which the Poisson equation for the linear displacement 
potential is evaluated at the initial Lagrangian position; see the accompanying paper [19]. 



22, 27], 



V q -^ 1 \q,T) = -D(T)6^(q,To) = -D(r)6o(q), 



(2- 



where D{t) is the linear growth function normalised to unity at z = 0, and the solution is 
purely longitudinal. Working with the solution (2.9), it is then easy to build iteratively from 
equations (2.7) and (2.8) higher order solutions. Expressed in Fourier space and keeping only 
the fastest growing mode, the resulting nth order displacement field Sb( n \p) = J-~[\I/( n )(g)](p) 
is then 
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(2vr) 3n 
x S (n) (pi, . . . ,p n ) <5 (pi) • • • 5 (Pn) , 



S (2) (pi,p 2 ) 

S (3) (pi,P2,Ps) 

S' (4) (Pl,P2,P3,P4) 



Pi 

2 ' 
Pi 



3^12 ^ 2 



.(») 



7 



,2„2 ' 



2 2 2 
Pl2 PlPi 

Pl23 

P123P1P2P3 
P1234 

Pi234PipIp1p4 



10 



3 3Q 21 36 



1 

+ - 



U> 



Ivi 



£> 



2 2 2 2' 

Pl23PlP2P3 



1 



+ 



51 U) 13 («) 14 ( a ) 20 ( s ) 1 ( a ) 

539 4a 154 4b 33 4c + 33 4d + 11 4e 

-J. s ) 



1 [1 (s) _ ]>_ (s) J_ (5) 

p? 2 34p!pIpIp! ie^ 21 w * + u w « 



(2.10) 

(2.11) 
(2.12) 
(2.13) 

2.14) 



.(») 



»/ 



and wi 



Here, the symmetrised kernels ftn' = «n (pi, • • • ,Pn) and u>n' = ujn'{pi, ■ ■ ■ ,p n ) represent 
respectively the longitudinal and the transverse component of if>( n >. Their exact forms can 
be found in appendix A.l. We emphasise again that while the displacement field can be split 
into purely longitudinal and purely transverse components, longitudinal-transverse mixing 
does occur generically when \I/ is used to reconstruct the density contrast via equation (2.3). 
As we shall see later, for our particular application, there is a mixing at fourth order between 
u^J and S^ ' (see also in appendix F). In more general cases, however, longitudinal-transverse 
mixing happens already at third order [29]. 

2.2 Bispectrum in the Lagrangian framework 

Starting with the Fourier transform of the density contrast [30], 



S(k) 



d 3 xe ikx 



S(x) 



d 3 qe 



\kq ( e ik*(q) 



1 



(2.15) 



where we have used equation (2.3) for the last equality, the bispectrum B{k\, k%, £3) can be 
defined as 

5(fci)5(fc 2 H(fc 3 )) =(2vr) 3 4 3) (fci + fc2 + fc 3 ) B(h,k 2 ,k 3 ). (2.16) 

/ c 

Here, hi = \h%\, the subscript "c" denotes the connected piece, and the Dirac delta follows from 
statistical homogeneity, so that the three wavevectors fei,2,3 always form a closed triangle. 
The dependence of B(ki, ^2,^3) on only the magnitude of the wavevectors is a consequence 



of statistical isotropy. Physically, this means that the bispectrum depends only on the shape 
and the size of the triangle formed by fei,2,3, not on the triangle's orientation. Keeping this 
in mind, we can rewrite the bispectrum (2.16) using the last equality in equation (2.15) as 

B(h, k 2 , fa) = /d 3 A 21 /" d 3 A3ie ifc 2 -A 21 +ifc3-A 31 

X /7 e ife 2-(*(q2)-*(<7i))+ifc 3 -(*(<?3)-*(qi))\ _ A ; (2.17) 



where Ajj = Qj — qj. Note that in deriving equation (2.17), we have used (expjifc • *&}) = 1, 
because of (S(k)) = 0. 

To evaluate the ensemble average in equation (2.17), we use the cumulant expansion 
theorem [15, 31], 

f °° -N ) 

(e iX ) = exp £ — (X N ) C , (2.18) 

Iat=i ' J 

where, in our case, X = k 2 ■ [^(q 2 ) - *(<?i)] + &3 • [*(93) - *(9l)]> an d (X N ) C denotes the 
A^th cumulant. Two options are available to us at this stage: 

1. A brute-force Taylor expansion of 8 in equation (2.15) up to the desired order in \I/, 

which yields 5 « {ft > + 5^ 2 ' + . . . for equation (2.16). 

2. We generalise the resummation scheme of [15] to the bispectrum calculation, and reor- 
ganise the perturbative series before performing a Taylor expansion up to the desired 
order in Vl/ '. 

As we shall see in section 3.4, option 1 in fact leads to the same outcome as standard 
perturbation theory (SPT), up to the same order in 5 (this is true at least up to fourth 
order). 

To implement option 2, we first rewrite equation (2.18) as 



E^<*">. -«p £>» «p E^E 

N=l ' ) KN=1 ' J |JV=1 ' j=l 



(2.19) 



where we have used the binomial theorem, and introduced the quantities 

D N ^(a N ) r + (b N ) r + (c N ) r , (2.20) 



E N = {W~*) c + J2 (l) {a k V- k c N -i) , (2.21) 

fc=o ^ ' c 



with the shorthand notation a = k\ ■ \l/(gi), b = k 2 ■ ^(92)1 and c = k^ ■ ^(q^). Using 
equation (2.19) and the fact that the Djy term is independent of Ay, the bispectrum expres- 
sion (2.17) now becomes 

B(h, k 2 , fa) = exp J f; ^D N (0) 1 j d 3 A 21 J d 3 A3ie ifc 2 -A 21 +i fc 3-A 31 

-{e5£(;)- 



^ - 1 , (2.22) 



where, because of statistical homogeneity, we can choose to evaluate Djy at any single point 
q in space, e.g., at q = 0. For large separations Ay, one expects the cumulants in Djy to 
have a larger contribution to B(k\, &2, ^3) compared with those in En [15]. Therefore, we 
expand perturbatively only the cumulants appearing in En, while the exponential prefactor 
is kept as it is, with the Dn terms in the exponent evaluated according to the Zel'dovich 
approximation [13, 22]. This is thus the "resummation" scheme of [15], and equation (2.22) 
is the starting expression for the so-called resummed bispectrum. 

3 Results I: Clustering in real space 

Working with the expressions (2.17) and (2.22), we now evaluate the bispectrum and the 
resummed bispectrum respectively up to one-loop order. For readers wishing to skip the 
technical details, the main result of this section is equation (3.15). A numerical comparison 
of our results with an "exact" bispectrum extracted from TV-body simulations can be found 
in section 3.3, while the equivalence between the LPT and the SPT bispectra is discussed in 
section 3.4. 

3.1 N-point correlators of the displacement field 

The terms Dn and En in equation (2.22) contain cumulants that are functions of the displace- 
ment fields *&(q) (see equations (2.20) and (2.21)). Since the perturbative solutions of *&(q) 
have been conveniently expressed in terms of their Fourier transforms in equations (2.10) 
and (2.11), we rewrite the cumulants similarly in Fourier space [15]: 

%M ■ ■ ■ *i N (PN)) c = (2vr) 3 4 3) (Pi-iv) C h ... iN ( Pl , ..., PN ). (3.1) 

Here, Cj 1 ...j JV (pi, . . . ,Pn) is an iV-point correlator, and its perturbative version is defined via 

«i ai) (Pi) • • • *t"W)) c = (2^) 3 4 3) (pi...iv) 4 a ™V, ...,Pn). (3-2) 

Note that, unlike the iV-point functions of the density contrast, these iV-point correlators 
are in general complex quantities. 

Assuming Gaussian initial conditions, the correlators relevant to our problem are: 

C«(p) = 4">(p) + Cf>(p)+Cf '"(p), 
C iJ ,( Pl , P2 , P3 ) = 4 211>,+C r+4i 411)) +»nC«r», (3.3) 

r (r, « « « ^ - ^{(1122)} ^{(3111)} 

where we have employed a shorthand notation for the 2-point correlator, Cij(p) = Cy(p, —p), 

and omitted writing out explicitly the dependences of C^, N = C^, N (pi, . • ■ ,Pn)- The 
curly brackets in the superscripts denote summation over all possible permutations, e.g., 

C ijk = C ijk + C ijk + C ijk ' wnile tne notation i ffi nCy fe = iC ijk + iiQjk indicates two 
distinct contributions to correlators of the type C^ k . In case the initial conditions contain 
PNG, the correlators 

C ij (p)=cf l)} (p), 



$ u ) + n c { '-' [i)] ../<i (22J <! , . ...H'- m )! 

C ijkl (pi,p 2 ,P3,P4) = C\]H l) + C ijk! 



C ijk (PuP2,Ps) = C^> + n C}r n + ieii^^ + leii^ 1 " , (3.4) 



1111) . ^{(2111)} 



must be considered in addition. 

We can now evaluate equations (3.3) and (3.4) using the perturbative solutions of ^ 
given in equations (2.10) and (2.11), and express each correlator in terms of the linear power 
spectrum P^^k, z) = D 2 (z)Po(k), where D{z) is the linear growth function at redshift z, 
normalised such that D(z = 0) = 1, and Po(k) is the present-day linear power spectrum. In 
the case of non-Gaussian initial conditions, the additional correlators (3.4) will also depend on 
the present-day linear bispectrum Bq and trispectrum To. The exact forms of all correlators 
in equations (3.3) and (3.4) can be found in appendix B, while in appendix C we give the 
relations of Bq and To to their primordial counterparts B§ and Tj> predicted by inflation. 

3.2 The LPT bispectrum and the resummed bispectrum 

We now use the perturbative correlators from section 3.1 to evaluate the (unresummed) 
bispectrum (2.17) up to one-loop order. We begin by separating the perturbative expansion 
of B into a tree-level part B^' and a one-loop contribution B^ 1 ' , i.e., 

B = B {0 '> +B^ + ... , (3.5) 

which are themselves split respectively into a Gaussian and a non-Gaussian contribution, 

2?<°> = B 211 + B , BW = £g) uss . an + 5 W Q , (3.6) 

where B 2 \\ is the tree-level bispectrum with Gaussian initial conditions, and Bq is present 
only in the case of non-zero PNG. The exact form of Bq is reported in appendix C. Here, we 
explicitly derive -B211 hi order to demonstrate how the formalism works. 

The bracketed term ((• • • ) c — 1) in the expression (2.17) can be Taylor-expanded using 
the cumulant expansion theorem (2.18) as 

(e iX ) c ~ 1 = 4<* 3 >c + ^(X% + \{X 2 ) C {X 2 ) C 

+ ^(X 2 ) C (X 3 ) C - ±(X 2 ) C {X 2 ) C (X 2 ) C + 0(X 7 ) , (3.7) 

where X = k 2 • [^(92) — *(<7i)] + &3 " [^(93) — *(Qi)]- Considering only Gaussian initial 
conditions, the lowest order contribution is 



.(11), 



B211 = -ik u k 2j k 3l Cgf ^'(fci.fca,**) + k 1 ik 2j C^ , (k 1 ) k 2l k 3m C^'(k 3 ) 

+ k u k2jC\] 1 \k 2 ) k u k 3m c\^\k 3 ) + fciifeyCg^fci) k 2l k 3m C^\k 2 ) , (3.8) 



,, v-2) hifomC^ (k 3 ) + k u k 3j c\] 
which, using the definitions of the correlators, can be evaluated to give 



B 211 = 2P L (k 1 )P L (k 2 )F 2 {s \k 1 , k 2 ) + two permutations , (3.9) 



with the kernel 



„{s) n ,, 5 lfci -k 2 (fa k 2 \ 2(ki-k 2 f 

^ u ^^) = 7 + 2^7U + w + 7n^r- (3 - 10) 

But equations (3.9) and (3.10) simply reproduce the results of SPT. Thus, at tree level, the 
LPT and SPT bispectra are identical. 



Similar expansions up to higher orders in \l/ lead to the one-loop expressions 

-^Gaussian = #41101230222 + #112201113 + #xx0yy + #110211 — #11011011 , , . 

(1) (3.11J 

-^PNG = #11201220113 + #111101112 + #11012 + #110111 , 

where explicit forms of the constituent terms are given in appendix D. The © symbol denotes a 
grouping of several terms with a common origin in the Taylor expansion, e.g., the constituents 

of #4ii©i23©222 = #411 + I©II#123 + #222, with i e n#i23 = I#123 + II #123, all originate from 
the {X 3 ) c term in equation (3.7). The ® symbol indicates a contribution from the "product 
terms" in the Taylor expansion, e.g., #ii<g>n®ii arises from (X 2 ) C (X 2 } C (X' 2 } C , and consists 
of products of three C- 's. 

In order to compare the LPT one- loop expressions with their SPT counterpart, we apply 
a diagrammatic technique which allows us to regroup the LPT contributions in terms of the 
diagrams they produce. In our notation, this means rearranging the LPT contributions into 
the groupings 

^Gaussian = S 411 + #123 + #222 , _ _ 

-^PNG = ^112 + #122 + #113 • 

Since SPT produces the same classes of diagrams, our LPT results can be compared with 
standard SPT results on a diagram-to-diagram basis. The various classes of SPT and LPT 
diagrams up to one loop are shown in figure 1. In appendix E we demonstrate how to 
construct diagrams in LPT, and report the regrouped terms #. Suffice to say, all regrouped 
terms in LPT agree with their SPT counterparts: 

B (o) _ B (o) B (i) _ B W ( o 13) 

-°LPT ~~ -°SPT ' ^LPT ~~ -°SPT ■ {O.iO) 

Since the expressions are identical, henceforth we shall omit the subscripts "LPT" and "SPT" . 
To compute the resummed bispectrum, we expand only those terms appearing inside the 
A21 and A31 integrals in equation (2.22), while leaving the exponential prefactor untouched. 
Expanding up to one-loop order, we find the resummed bispectrum 

B(ki, k 2 , k 3 ) = exp I -- / ——3 [kukij + k 2 ik 2 j + k 3i k 3j ] C^ 1 ' (p) 
B (0) (fci, h, k 3 ) + # (1) (^i, k 2 , k 3 ) + lB^(h, k 2 , k 3 ) 



2 
d 3 p 



7^3 [hikim + k 2[ k 2m + k 3l k 3m ] C\^\p) 



(3.14) 



Given C\- (p) = (piPj/p 4 )PL(p) from equation (B.l), the angular integration of the first and 



the last terms can be easily performed, thereby leading to 



(kt + m + ki) r 

B(h, k 2 ,k 3 ) = exp I - V ± ^ 3J J d P P L (p) 



#(°) +#« + ^^fj^ B^JdpPUp) 



(3.15) 






a.) b.) c.) d.) 







e.) f.) g.) h.) 






i.) j.) k.) 



Figure 1. Diagrams in SPT and LPT. The symbols 0, A and ■ denote respectively a linear power 
spectrum, linear bispectrum, and linear trispectrum. For Gaussian initial conditions, the contributing 
diagrams are (a) -B222, (b) 1-B321, (c) 11-B321, (d) B411, and (e) B 2 n- For non-Gaussian initial condi- 
tions, the additional contributions are (f) uB 2 n, (g) 1-B113, ( n ) 11-8113, (i) 1S122, (j) 11-6122, and (k) 
Bm = Bq. As discussed in the main text, the sum of all LPT contributions to one specific diagram 
leads to analytic agreement with SPT. 

with B^ ' and B^- 1 ' given in equation (3.6). Equation (3.15) is the main result of this paper. 
The resummed power spectrum was calculated in [15]. Comparing it with our resummed 
bispectrum, one immediately recognises an overall similarity, especially in the form of an ex- 
ponential suppression prefactor. Analogously, as we show in section 3.3, while equation (3.15) 
generally constitutes a better approximation of the bispectrum in the weakly non-linear 
regime, the highly non-linear regime is dominated by the unphysical damping factor. 

3.3 Comparison with N-body results 

In this section we compare our resummed bispectrum (3.15) with bispectra extracted from 
A-body simulations. We use the A"-body results of Sefusatti et al. 2010 [32], read off their fig- 
ures 1, 3 and 5 using the plot digitiser EasyNData [33]. The simulations have been performed 
in a box of side length 1600 /i -1 Mpc spanned by a 1024 3 grid, for a ACDM cosmology with 
parameters /i = 0.7, S7 m = 0.279, f)(,= 0.0462, n s = 0.96, and a fluctuation amplitude <ts = 0.81, 
with initial conditions generated at redshift Zi = 99 using the Zel'dovich approximation. For 
the numerical evaluation of equation (3.15), we have written a C++ code, wherein the integrals 
are evaluated using the deterministic integration routine CUHRE from the CUBA library [34]. 
The code takes as an input the linear power spectrum calculated with CAMB [35] , and linearly 
interpolates it for the purpose of the loop integration. 

Gaussian initial conditions. Six separate comparisons are shown in figure 2 for the case 
of Gaussian initial conditions, corresponding to bispectra of the equilateral configuration 
{k\ = /c2 = A/j = k) and a squeezed configuration {k\ = ki = k, k% = Ak = 0.012 h/Mpc), each 
at redshifts z = 0,1,2. Each bispectrum displayed has been normalised to the smoothed, 
no- wiggle tree- level bispectrum B nw computed using the transfer functions of [36]. In each 
panel, the green dotted line denotes the SPT result without resummation, i.e., B(°> + B^- 1 ' , 
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Figure 2. Comparison of the one-loop corrected matter bispectra from SPT (green dotted line) and 
resummed LPT (RLPT; solid red line), with the iV-body results of [32] (data points), at z = (first 
row), z = 1 (second row), and z = 2 (third row). Gaussian initial conditions have been assumed. 
Left: Bispectrum of the equilateral configuration. Right: Bispectrum of one squeezed configuration, 
with A/c = 0.012ft,/Mpc. The red arrows indicate the regions of validity for the RLPT results at 
each redshift. For reference, we show also the tree-level bispectrum (black dot-dot-space line). All 
bispectra have been normalised to the no-wiggle tree-level bispectrum. 
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Figure 3. Comparison of the one-loop corrected matter bispectra from SPT (green dotted line) 
and resummed LPT (RLPT; solid red line), with the ./V-body results of [32] (data points) at z = 0, 
assuming PNG of the local kind and /nl = 100. Left: Bispectrum of the equilateral configuration. 
Right: Bispectrum of one squeezed configuration, with Afe = 0.012/i/Mpc. In each case, the ratio 
B(/nl = 100)/B(/ nl = 0) is shown. 

while the resummed bispectrum (RLPT) of equation (3.15) is represented by the red solid 
line. 

A general deviation of RLPT from the ./V-body results is expected in the highly non- 
linear regime, because the exponential prefactor in equation (3.15) has been evaluated only 
within the Zel'dovich approximation. Therefore, in order to compare our RLPT results with 
./V-body data, we must first define a cut-off scale k e s [15], 



k cS = a 



1 

12^ 



dpP L (p) 



(3.16) 



beyond which the exponential damping factor becomes too efficient for our RLPT results to 
remain physical. Note that k e g is time-dependent, and scales with the linear growth factor 
D(z) as \/D{z) per definition. The parameter a is a fudge factor that must be adjusted to 
the TV-body data. For a = 1/3, we find good agreement between RLPT and iV-body data in 
the k < fc e ff region of validity (indicated by the red arrows in figure 2). The high redshifts 
results are especially encouraging: at z = 2, RLPT remains compatible with iV-body data up 
to k ~ 0.2 h/Mpc, and appears to provide a better approximation to the iV-body bispectrum 
than does SPT. 

Unfortunately, an error estimation of our approximation with respect to the iV-body 
data is in general not possible, because the simulation errors are often very large. This is 
especially so in the case of the squeezed triangle, where the short side Ak = 0.012 /i/Mpc 
is merely a factor of three larger than the fundamental wavenumber of the simulation box 
(~ 0.004 /i/Mpc). Sampling errors are expected to be large in this instance. 

Non-gaussian initial conditions. Figure 3 shows the case of non-zero PNG of the local 
type, characterised by /nl = 100 (see appendix C for the definition), and its effect on the 
matter bispectrum at redshift z = 0. We consider again the equilateral configuration and the 
same squeezed configuration as above, computed from one-loop SPT, resummed LPT, and, 
for reference, tree-level SPT. 

Following [32], we plot the ratio -B(/nl = 100)/.B(/nl = 0), where /nl = corresponds to 
Gaussian initial conditions. Results from one-loop SPT are represented by the green dotted 
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line, while the red solid line shows our RLPT calculations. The reference tree-level results — 
B211 for Gaussian initial conditions, and B211 + Bq for the non-Gaussian case — are denoted 
by the black dot-dot-space line. For convenience we use the non-Gaussian contributions at 
tree-level and one-loop from figure 1 of [32]. 

Since the PNG considered here is of the local type, we expect its contribution to the 
present-day matter bispectrum to be peaked in the squeezed configuration. Indeed, the 
contribution due to PNG in the mildly non-linear regime (k ~ 0.04-0.1 /i/Mpc) is at the 
11% level in the squeezed case, while for the equilateral configuration the contribution in 
the same k range is about 3%. The SPT result shows a slight increase in the importance of 
PNG at higher /c-values, while the RLPT result suggests a constant behaviour. Comparing 
with TV-body data, we find that the iV-body results tend to overshoot both the RLPT and 
SPT results. It is not clear to us at this stage whether this is a problem of the semi-analytic 
methods, or of the iV-body simulations. We note however that one of the main hurdles facing 
simulations with PNG is the accurate implementation of non- Gaussian initial conditions, and 
much research in this direction is ongoing (e.g., [37-39]). We conjecture that the generation 
of initial conditions may contribute at least partly to the discrepant results. 

3.4 Exact relationship between SPT and LPT 

We have demonstrated in section 3.2 and the associated appendices that the one-loop bispec- 
tra computed from SPT and LPT are identical. Similarly, the equivalence between the SPT 
and the LPT one-loop power spectra was previously shown in [15]. Here, we summarise a 
result from our accompanying paper [19], which shows that LPT and SPT in general return 
the same results at least up to fourth order in the density contrast 5. 

Our starting point is equation (2.15). To prove the equivalence of LPT and SPT, we 
expand the LHS of equation (2.15) as 

S = sW + 6^+m+6^ + ..., (3.17) 

and similarly the exponential on the RHS up to the same order in the displacement field \I/. 
Using explicitly the fastest growing solutions for the nth order displacement fields given in 
equation (2.10) in the initial position limit, and summing up all contributions with n powers 
of 5q, we arrive at 

**>(*,.) = ^/ d3p ^ ) s d > (2^4 3) (pi--.-fe) 

x X^(k; Pl , . . . , p„) S„( Pl ) • • • 6 {Pn) , (3.18) 

(s) 

where J„ is a symmetrised kernel at nth order. Explicit forms can be found in appendix F 
for n < 4. The Dirac delta in equation (3.18) fixes k = p\... n . A careful examination of the 

(s) 

kernels X n then reveals 

X n s \k; Pl ,..., Pn ) =4 s) (pir--P«), (3-19) 

k=p 1 ... n 



where Fn is simply the nth order SPT kernel for the nth order density contrast [12], shown 

(s) (s) 

in a symmetrised form in appendix A. 2 for n < 4. The equivalence between Xn and F n 
has been checked algebraically with a Mathematica code up to n = 4 [19], which is available 
upon request. 3 

E-mail to rampf@physik.rwth-aachen.de. 
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4 Results II: Clustering in redshift space 

The resummation technique used in this work can be easily generalised to the calculation of 
the matter bispectrum in redshift space [40] . We briefly outline the computational procedure 
in the following. The final expression for the resummed redshift-space bispectrum up to one 
loop, assuming Gaussian initial conditions, can be found in equation (4.16). 

4.1 Density contrast and distorted displacement field in redshift space 

In the so-called plane parallel limit, the comoving distance s in redshift space is defined 

as [15, 41] 

s = x + — -z = x-u z z, (4.1) 

where x is the comoving distance in real space defined in equation (2.1), z is the line of sight 
(assuming the observer is fixed on the comoving grid), and v = dx/dr is the peculiar velocity 
of the fluid element. Rewriting equation (4.1) as 

s(q,T)=q + * s (q,T), (4.2) 

it is easy to see that the distorted displacement field reads 

* S = *H i — z . (4.3) 

7x 

Mass conservation in redshift and Eulerian space implies [42] 

p(l + 5 s {s))d 3 s = p{l + 6{x))d 3 x = pd 3 q, d 3 s = J s d 3 q, (4.4) 

which leads immediately to 5 s (s) = 1/J S — 1, and its Fourier transform 

5 s (k) = fd 3 qe ikq ( e ik -* s ^ - l) . (4.5) 

The redshift-space power spectrum and bispectrum can then be defined as 

[S s {k)6 s (k')\ =(2ir) 3 5 ( n\k + k')P s (k), (4.6) 

^ s (fei) 5 s (k 2 ) 5 s (fc 3 )) c = (2vr) 3 4 3) (fci + k 2 + fc 3 ) B s (fci, fc 2 , fe 3 ) , (4.7) 

on which we can now apply a perturbative procedure analogous to that outlined in section 2.2 
to obtain tractable results. 

4.2 The resummed bispectrum in redshift space 

To compute the perturbative correlators in redshift space, we first note that the perturbative 
kernels S^ n ' given in equation (2.11) are time-independent; the whole time evolution of the 
displacement field \l/( n ) is embedded in a separable linear growth function D as \I/( n ) oc D n . 
Therefore, the time derivative of the nth order displacement field can be written as 

^- = nUf* (n) , (4.8) 
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where / = dlnD/dlna is the logarithmic growth rate. From equation (4.3) it follows that 

or, in index notation, 

(4-10) 



**("> = (5 tJ + nfZiZj)^ = i^ n) 



where, again, summation over repeated indices is implied. 

The tensor R^ deforms the initially isotropic displacement field \I/' n ), thereby ren- 
dering its redshift-space counterpart <|> s ( n ) in general direction-dependent. It follows that 
the perturbative correlators constructed from the displacement fields must also be similarly 
distorted, i.e., 



a 



s(a±---a N ) , 



? (ai) n(a N ) ^(cLf-aN) , 



(4.11) 



and computing the redshift-space bispectrum becomes simply a matter of replacing all oc- 
currences of Cj-"^ aN with C\ ' „^ . For example, the tree-level redshift-space bispectrum 
assuming Gaussian initial conditions is given by 

B s 2X1 (k 1 ,k 2 ,k 3 )= (4.12) 

-ik u k 2j k 3l Cg (211)} (fc!, k 2 , k 3 ) + kuk^C^ikx) k 2l k 3m C^ l \k 3 ) 



-."(n)/ 



lS (n) 



"(ii). 



5(H)/ 



+ k u k 2j C* xx, (k 2 ) kuhmC^'iks) + kuksjC^'ikx) k 2l k 3m C^'(k 2 ) , 



cf the real-space expression (3.8). Using the relation (4.11), we obtain, after some reshuffling, 
the result 



B s 2n (k h k 2 ,k 3 ) = 2C(k 1 ,k 2 ,k 3 )R(k 2 )R(k 3 )P L (k 2 )P L (k 3 ) 

+ two cyclic perms. , 



(4.13) 



where 



C(fti, fc 2 , fe 3 ) = F 2 (s) (fe 2 , fe 3 ) + f {kl k2 z) G 2 s \k 2 , k 3 ) 



^f(kx-z) 



{k 2 -z) f 1 + f (k 3 -z) 



k 2 



k?, 



+ 



(k 



3 • Z) 

k 3 



R(ki) = 1 + / 



{kj ■ z) 



1 + / 



'M ■ z) 

k~± 



(4.14) 



(4.15) 



with the symmetrised second-order SPT kernels F 2 and G 2 (see appendix A. 2). The 
redshift-space kernels R(ki) and C(k\,k 2 ,k 3 ) agree exactly with the results from SPT; in 
the notation of [42], they are called Z\ and Z 2 , respectively, with 6 = 1 and b 2 = (i.e., 
neglecting local galaxy biasing). 

The one-loop correction to the redshift-space bispectrum can be computed in a similar 
manner. Generalising the procedure to the resummation case, we find the following expression 
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for the resummed redshift-space bispectrum up to one loop: 



B s (ki,k 2 ,k 3 ) = exp < 



£L kf + f(f + 2)(k t -z) 



12ir 2 

>3 ,,2 



dpP L (p) 



(4.16) 



s «(o) + B *(i) + 



n U=i[kt + f(f + 2)(ki-zy 



12vr 2 



£ s ( ) / dpP L (p) 



Note that, unlike the case of the real-space bispectrum, redshift-space distortions induce for 
B s a dependence on the direction of the constituent wavevectors. We leave the numerical 
evaluation of equation (4.16) for a future project. 

5 Conclusions 

In this paper, we have computed the matter bispectrum up to one loop using Lagrangian 
perturbation theory (LPT). We find an exact agreement between the bispectrum obtained 
this way and its counterpart from standard Eulerian perturbation theory (SPT), as long 
as the LPT calculation is restricted to the initial position limit [19]). More generally, in 
this limit, we find that both LPT and SPT predict the same density contrast up to fourth 
order. The analytic agreement between these theories up to third order can be shown by 
hand, while the equivalence of the fourth-order solutions can only be established algebraically 
using a Mathematica code (which can be obtained upon request). Note that the use of the 
initial position limit is, to our knowledge, implicit in all LPT calculations in the current 
literature performed in Fourier space (e.g., [13, 15, 16]). Relaxing this assumption might 
yield better approximations of the matter TV-point spectra in the weakly non- linear regime, 
and should be a promising avenue to study in the future. 

We have also generalised a resummation technique, first developed in [15] to compute a 
resummed matter power spectrum, to the matter bispectrum calculation both in real space 
and in redshift space. In the mildly non-linear regime where we expect the method to be valid, 
we find good agreement between the "resummed" real-space one-loop bispectrum and the 
"exact" bispectrum extracted from the iV-body simulations of [32] assuming Gaussian initial 
conditions. Our results for bispectra of the equilateral configuration at high redshifts (z = 2) 
are especially encouraging: within the region of validity, the resummed bispectrum (RLPT) 
generally performs better than its one-loop SPT counterpart in terms of their comparison 
with the iV-body data. The performance of the resummation scheme is more difficult to 
judge for bispectra of a squeezed configuration. Here, for a squeezed triangle with a short 
side Ak = 0.012/i/Mpc, we find that the TV-body bispectra of [32] at all considered redshifts 
have more power in the k region of interest, compared with all semi-analytic approximations 
studied in this work (i.e., tree- level, one- loop SPT, and resummed one- loop). However, 
we note also that the Ak value used in the comparison is rather close to the fundamental 
wavenumber of the simulation box (~ 0.004 /i/Mpc). This suggests that the iV-body bispectra 
in the squeezed limit could very well be afflicted by large sampling errors. In the case of 
non-Gaussian initial conditions of the local variety (with /nl = 100), we observe a slight 
overshoot of the iV-body results of [32] compared to RLPT (and often to SPT as well). 
Further studies will be required to understand this discrepancy, but we note for the moment 
that the crucial issue of generating non-Gaussian initial conditions in iV-body simulations is 
highly non-trivial, and research in this direction is ongoing (e.g., [37-39]). 
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The consideration of biasing effects in the galaxy clustering statistics (e.g., [43-45]) is 
beyond the scope of this work. However, a simple replacement of / — >• /? = f/b [15] in 
our resummed redshift-space bispectrum (4.16) immediately generalises the expression to 
describe also a simple linear biasing scheme with bias parameter b. Extending this scenario 
to non-local biasing schemes is less straightforward, and will be addressed in a future work. 
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A Perturbative kernels 

A.l LPT 

The unsymmetrised kernels n n and u n in equations (2.11) to (2.14) are given by 

K2(Pi,P2) = g b?P2 - (Pi -Pa) 2 ] , (A.l) 

K3a(Pl,P2,P3) = -g [P1P2P3 ~ (Pi ■ Plfpl ~ (Pi ■ Pzfpl ~ (P2 ■ Pzfvl 

+2( Pl ■ p 2 )(pi ■ p 3 )(p 2 ■ Pa)] , (A.2) 

I \ t s K 2(Pl,P23) ff . Q x 

«36(Pl,P2,P3) = -K 2 {P2,P3) 2 ' A - 3 

P23 

I \ I \ ( \ K 2(Pl3,P24) , A .,. 

«4a(Pl,P2,P3,P4) = K2(Pl,P3)K2(P2,P4) 2"^ ' ( A - 4 ) 

^13^24 

K 2(P3)P4) 

«46(Pl,P2,P3,P4) = -6K3o(Pl,P2,P34) 2 ' ( A - 5 ) 

P34 

/ s / s«2(Pl,P234) , A c s 

K4c(Pl,P2,P3,P4) = -«3a(P2,P3,P4) 2 , (A.6) 

^234 

/ s / ^ 2 (P1,P234) fK _v 

«4<i(Pl,P2,P3,P4) = -K 3 &(P2,P3,P4) 2 ' ( A - 7 ) 

P234 

K4e(Pl,P2,P3,P4) = ~ [(Pi • P2)(P234 • P34) ~ (Pi • P34)(P234 • P2)] 

X (Pi -P234) -| f~ «2(P3,P4) , (A.8) 

P234%4 



and 



W3c(pi,P2,P3) = [Pi (P123 • P23) - P23 (P123 • Pi)] (Pi • P23) — |p~^~ ' ( A - 9 ) 

^23 
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W 4 /(Pl,P2,P3,P4) = [Pi (Pl234 • P234J ~ P234 (Pl234 ■ Pi)] (Pi ■ P234J 2 

P234 



<^4g(Pl,P2,P3,P4) = [Pi (P1234 • P234) ~ P234 (Pl234 " Pi)] (Pi " P234; 



^3fe(P2,P3,P4) 



^234 



<^4h(Pl,P2,P3,P4) = [Pi (P1234 • P234) ~ P234 (Pl234 "Pi)] 

The symmetrisation procedure is then straightforward [46]: 



» 



Pl-^ 3c (P2,P3,P4) 



^234 



(A.10) 



.,(-) 



(Pi, • • • ,Pn) = — Y^ K n(Pi(l),- ■ ■ ,Pi(n)) 



ieS n 



U 



w 



(Pi, • • -,Pn) = —^2 W n(Pi(l)' • • • 'Pi(n)) • 



(A.ll) 



ies„ 



Note that the kernels u satisfy the condition 

Pi2-n ■u) n (p 1 ,...,p n ) = 0, (A.12) 

indicating transverseness. 

A.2 SPT 

The functions F n and G n arising from a perturbative solution of the Eulerian fluid equa- 
tions (2.4) to (2.6) in an Einstein-de Sitter universe can be obtained using a set of recursion 
relations [12]. Explicit expressions up to n = 4, in their unsymmetrised forms, can be found 
in, e.g., [46]. Here, we give the symmetrised versions of these expressions. 



For n = 2 



.,„,. .. blk l -k 2 (k 1 k 2 \2{k 1 -k2) 2 



i{s), 



G 2 s '(kt,k 2 ) = --\- 

For n = 3: 

F( s \k u k 2 ,k 3 ) 



3 , lfci-fca ^1 fo\ , 4(fci-fc 2 ) 2 



7 ' 2 A:iA; 2 \k 2 ' hj ' 7 fcfjfcf 



1 l fc 12 3-fc 3 ( S ) 

3 I 18 JfeS 2 V ; 



+ 



7 fcl23 • ^12 1,2 fcl2 • &3 

9 ^1^3 



18 fc 2 2 



G 2 s) (fei,fc 2 ) 



(A.13) 

(A.14) 



2 cyclic permutations > , (A. 15) 



G^(k lt k2,k 3 ) 



1 fcl23 -fc3 p (s). 



18 fc 2 2 l l! 2j 



1 fcl23 • fcl2 1 ,2 fclg • &3 

18 k 12 9 ^i 2 ^3 



G^(fei,fe 2 ) 



+ 2 cyclic permutations > . (A. 16) 
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For n = 4: 



,2 

'1 



_fel9.S4 • fcl 



396-F 4 w (fci,fc 2 ,fc 3 ,fc 4 ) = <J27 ^ - F 3 w (fc 2 ,fc 3 ,fc 4 ) 



+ 



97 fcl234^fc234 2 fcj • fc 234 

Z ' p " + DK 1234 ,21.2 

^234 ^1^234 



Gg (&2, ^3) ^4) + 3 cyclic permutations > 



+ < 18 -2 ^2(^1,^2)^2 f (&3, fc 4 ) + 5 permutations > 

I ^12 " J 

+ (^1234 ^2 ',t 34 G 2 S) ( fc i» fc 2)g 2 s) (fc 3 ,fc4) + 2permutationsl , (A.17) 

396-Gi a) (fci,fe 2 ,fe 3 ,fc 4 ) = |9 fcl2 ^ 2 ' fcl F 3 (a) (fe 2 ,fc 3 ,fc4) 

G 3 (fc 2 , fc 3 , fc 4 ) + 3 cyclic permutations 



n fcl23 4 • ^234 . ,2 ^1 ' fc 234 
y P I" Z4ft 1234^2p 

^234 ^1^234 

+ < 6 -2 G 2 (fci, fc 2 )F 2 (fc 3 , fc 4 ) + 5 permutations > 

+ { 16fc?234 ^ , 2 34 G' 2 S) ( fc i 1 fc 2)g 2 s) (fc3,fc4) + 2 permutations 1 . (A.18) 

I ^12^34 J 

The advantage of using the above expressions is that the symmetrised nth-order kernels have 
the property that, in common applications, partial sums of several wavevectors equate to zero. 
For instance, to compute the loop correction P13 to the power spectrum, one integrates over a 
kernel comprising the function F% f (q, —q, k). Plugging the arguments into equation (A. 15), 

we see immediately that the second term is proportional to G 2 (q, —q)/(q — q) 2 , where both 
the numerator and the denominator evaluate to zero. A careful treatment shows that this 
term in fact vanishes as e 2 • 1/e as e — > 0, which is not very straightforward to see from 
inspecting the corresponding expression in [46]. 

B Perturbative N-point correlators 

Here we give explicit expressions for the perturbative correlators (3.3) and (3.4), computed 
following the definition (3.2) of C^ 11 '-' and using the expressions (2.10) and (2.11) for <I>W. 
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B.l Gaussian initial conditions 

For Gaussian initial conditions, the correlators (3.3) are 

C% 1 \p)=S ( . 1 \p)sU(p)P L (p), (B.l) 

& 22 \p) = 2 J 0fSP(p',p- p>) Sf\p\p- p>) P L (\p - p'\) P L (p') , (B.2) 

C\f\p) = 3Si 1 \p)P L (p) J ^sf\p 1 p',-p')P L (p') = cf\p), (B.3) 

C y ijk '{Pl,P2,P3) = C] ki '(p2,p3,Pl) = C y kij '(P3,Pl,P2) 

= -2iS {2 \p 2 ,p 3 ) sV(p 2 ) P L ( P2 ) S^ips) P L ( P3 ) , (B.4) 

C^ 2 \puP2,P3) = 8iJ0^sl 2 \p 1+P ', -p')sf\- Pl - p',p'- P3 ) 

x PL(\P3-p'\)PUp)S { k 2 \p',p 3 -p')P L (\pi+p'\) , (B.5) 

C ijk (Pl,P2,P3) = C^ k ] i) (p2,P 3i Pl) = C£j 1 \p 3 ,p 1 ,p 2 ) 

= \2iSf\p 2 ) P L (p 2 ) 5«(p 3 ) PM J (§13^ 4) (-P2, -P3,P', -P') Pl(p') , (B.6) 
iC% 2 k 3 \pi,p 2 ,P3) = 6i^ (1) ( Pl ) P L ( P1 ) J J^ Sf\p\p 2 - p>) 

*PL{\P2-p'\)Sf\-p U p'-p 2 ,-p')P L (p'), (B.l) 

uCl- 2 k 3 \pi, P 2,P3) = -6i5| 1) (pi)P L (pi)5J 2) (PUP 3 )P L (P3) 

* j ^3 s k\p^p'l-p') p L(p')l (B - 8) 

IeIlCiJfc 3) (Pl'P2,P3) = lC^ k 3 \pi,P2,P3) +IlC[] k U \pi,P2,P3) 

= miC^jlf'(p 2 ,p3,Pi) = miC ( k ^ 2 \p 3 ,pi,P2) = miCl k 1 f\pi,p 3 ,P2) 

= miC k 21) {P3,P2,Pi) = miCj 2 l 3) (P2,Pi,P3) , (B.9) 

^(1122), x ^(1221)/ x ^(2211)/ s 

C; jU '(Pl,P2,P3,PA) = C y jkH '{p 2 ,P3,P4,Pl)=C kHj '(P3,PA,Pl,P2) 

/^(2112)/ x ^(2121)/ x ^(1212)/ x 

= C lijk (P4,Pl)P2,P3) = C likj m,Pl,P3,P2) = C ikjl (PUP3,P2,Pi) 

= -ASf\p 1 )P L {p x )sf\p 2 )P L {p 2 ) (B.10) 

x 's i k 2 \-p 1 ,p 13 )sl 2 \p 2 ,p 13 )P L ( Pl3 ) + ^ 2) (-P2,P23)sJ 2) (pi,p 23 )Pl(P23) 

^(1113)/ x ,-,(1131), x 

Cl jkl '{Pl,P2,P3,P4)= C y jkH '{p 2 ,p 3 ,P4,Pl) 

L1 )^ 

[P3,PA,Pl,P2 



C 'feS 11) (P3,P4,Pl,P2) = Cjjfc (P4,Pl,P2,P3) 



6 S^'( Pl ) Sy>(p 2 ) S ( k '( P3 ) Sl '( Pl ,p 2 , P3 ) P L ( P1 ) P L (p 2 ) P L (p 3 ) . (B.ll) 
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B.2 Non-Gaussian initial conditions 

As shown in equation (3.4), additional terms arise in case the initial conditions contain non- 
Gaussian contributions. 

c\f\p) = Cf\p) = J ^StPip) Sf\p\p- p') B (p,p>, \p -p'\), (B.12) 

c\]i 1 \puP2,Pz) = iSf\p 1 )sf\p 2 )s£\p ?J )B () (p 1 ,p 2 , P z), (B.13) 

nCjfc (Pl,P2)P3) = llCfki (P2,P3,Pl) = IlC2 J 11) (P3,Pl,P2) , 

= i5f 1) (pi)5J 1) (p 2 ) J 0^S^(p',p 3 -p')U Pl ,p 2 ,p',p 3 -p'), (B.14) 

iCyfc (Pl»P2>P3) = 3iB (Pl,P2,P3) 

x I ^^V) ^ (1) (P2) S£\p 3 ,p>, -p>) P L (p>) , (B.15) 
nCg 3) (pi,p 2) p 3 ) = -3iP L (pi) J0^sl 1 \p 1 )sf\p 2 ) 

XS^\p 1 ,p',p 2 -p')Bo(p2,p',\p2-p'\)-(pi<r+P2,i*+j) , (B.16) 

IeIlCiJfc 13) (Pl'P2,P3) = lC^ fc 13) (pi,p 2 ,P3) +IlC§ fc 13) (pi,p 2 ,P3) 

= ieiiCJii (P2,Ps,Pi) = i®nCj2j- (P3jPi.ps) . ( B - 17 ) 

iCg 2) (pi,p 2 , P3 ) = -2i sl 1 \ Pl )P L ( Pl )sf ( Pl , P3 ) 

X y ^^ 2) (P',P3 -P') B (p', P3 , \p'- P3 \) - (P3 ++ ft, j <-> fc) , (B.18) 

iiC^ 22) (pi,p 2 ,P3) = -AiS^\ Pl ) J ^Sf\p',p 2 -p') 

x ^ 2) (pi+j/,i>2-p')i'£(|P2-j/|)So(Pi,p',bi+y|) ) (B.19) 

^(122)/ n --.(122)/ v . ,-,(122), v 

I®nCyfc (Pl.P2,P3) = lQj fc (Pl,P2,P3)+IlC^ fe HPl,P2,P3J 

= ienC]S 1) (P2,P3,Pi) = ieiiCfc 2 j 12) (P3,Pi,P2) , (B.20) 

^ 11) (pi,p 2 ,p 3 ,P4) = Sf\ Pl ) Sf\p 2 ) s£\p 3 ) 5«(p 3 ) S z (1) (p 4 ) 

xT (pi,p2,P3,p 4 ), (B.21) 

^(1112)/ s .-,(1121), x 

C" ijkl (pi , p 2 , P3 , P4j = Cj M . (p 2 , P3 , P4 , Pi ) 

.-,(1211), x ^(2111)/ «. 

= C klij \P3,P4,Pl,P2) = C} ijk ; (P4,P1,P2,P3J 

= -2Sf\p 1 )P L { Pl ) sf\p 2 ) S£\p 3 ) sl 2 \p llP23 ) B (p 2 ,p 3 ,p 23 ) 
-2Sl 1 \p 1 )S^(p 2 )P L ( P2 )si 1 \p 3 )SP(p 2j p 13 )B (p 1 ,p 3 , Pl3 ) 
-2Sl 1) (p 1 )S^(p 2 )S^(p 3 )P L ( P3 )sl 2 \p 3 , Pl2 )B (p 1 ,p 2 , Pl2 ). (B.22) 
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C From primordial curvature perturbations to non-Gaussian initial con- 
ditions 

Primordial non-Gaussianity from inflation is usually quantified in terms of the iV-point statis- 
tics of the superhorizon metric perturbations <&(x), extrapolated to the epoch of matter 
domination. For our bispectrum considerations in section 3, up to one- loop order only the 
primordial bispectrum B$(ki, k 2 , k 3 ) and the primordial trispectrum T$(ki, k 2 , k 3 ,ki) play 
a role. To make use of these initial conditions, however, we must first link B§ and T$ 
respectively to the linear matter bispectrum Bq and trispectrum To at redshift zq via 

Bo(ki, k 2 , h) = M(ki,z ) M(k 2 , zq) M(k 3 , z ) B*(ki, k 2 , h) , 

r (fei, k 2 , fe 3 , fe 4 ) = Af (fci, zq) M{k 2 , z ) M(k 3 , z ) M(h, z ) T*(fci, k 2 , k 3 , fc 4 ) , (C.l) 

where 

M(M) = ?^^. (C.2) 

Here, the linear transfer function T(k) asymptotes to unity as k — > 0, while the linear growth 
function D(z) is defined to coincide with the scale factor a during matter domination. Note 
that this definition of D(z) differs from our previous definition in sections 2 and 3. We 
compute M(k,z) in this work using the publicly available Boltzmann code CAMB [35]. 

Assuming primordial non-Gaussianity of the local type, one can further parameterise 
the degree of PNG in terms of the parameter /nl> defined via [3] 

$(*) = ip(x) + /nl [p\x) - {ip 2 {x)j\ , (C.3) 

where <p(x) is a Gaussian field. The primordial bispectrum and trispectrum of <£ are then 
given by 

B$(ki, k 2 , k 3 ) = 2/ NL P$(fci) Pq,(k 2 ) + two perms. , 

r«(fci, k 2 , fc 3 , fc 4 ) = 4/^ P 9 (ki) P*(k 2 ) [P^(k 13 ] + P 9 (k 14 )] + 5 perms. , (C.4) 

to leading order in /nl- 

D One-loop expressions for the bispectrum 

D.l Gaussian initial conditions 

We split up the one- loop contribution to the matter bispectrum into five groups: 

-^Gaussian = -^41101230222 + -E>1122©1113 + -Bxx®y y + -Bll0211 ~ -Bll011011 , (D-l) 

where 

-641101239222 = -ikuk 2 jk 3k C>^ llel23 ® 222)} (fei, fe 2 , k 3 ) , (D-2) 

d 1 1 r 1 r f d 3 P /-rf (112201113)}/ , 1 j \ 

-D11220H13 = ^hikijk 2 kk 3 i / ——3 C>> u (p, fei - p, k 2 , k 3 ) 

. 1, , , , I d 3 p „{(1122©1113)}/ , , , n 

+ ^k 2i k 2j k lk k 3 i / -—3 C^ v fc , ;/ (p, fc 2 - P, fci, ^3) 

+ ^k 3i k 3j k lk k 2 i / -—3 C^ 122elU3)} (P> fc 3 - P, fci, ^2) , (D.3) 
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B 



xxftOyy 



kuk 2j k 2 ik 3m 


c\] 1 \k l )cl^ 22) \k,) + cf l ^ 2)} 


kuk 2 jk 1 ik :im 


'cf\k 2 )C^\k 3 ) + C f 1 ^ 


kuk 3 jk 2 ik 3m 


'cf\k 1 )c^ lm2) \k 2 )+cf lm2)] 



(k2)cH\k ?> ) 
(ki)C^\k 2 ) 



(D.4) 



#n®2ii =-- ^kiik2jk 3k c}jl n ' i (ki,k 2 ,k 3 ) 



2 
d 3 p 
(2^)3 

d 3 p 



(n)/ 



[kukim + k 2 ik 2 m + k 3 ik 3m ] C\ m '(p) 



+ i / jT—^kuk^j C\f-\p) ki k k 2[ k 3m c|^ 11)} (fci -p,p + k 2 , k 3 ) 
+ i / 77^3 fciiAty C^\p) ki k k 2L k 3m C|^ 11)} (fci - p, k 2 ,p + fc 3 ) 
+ i / jT^foifoj Clf\p) ki k k 2 [k 3m C|}^ 11)} (fci, k 2 - p,p + k 3 ) 



+ \k li k 2j c\f\k 1 ) 



+ ^k ll k 2j C\f\k 2 ) 



+ ^k lt k 3j ci; L) (k 1 ) 



+ ^k lt k 3j ci; L) (k 3 ) 



+ ^k 2i k 3j C\f\k 2 ) 
+ \k 2i k 3j clf\k 3 ) 



7^3 [k 3k k 3 ik 2m - k 2k k 2 ik 3m ] C^ 11)} (p, k 3 - p, 
777-T3 [k 3k k 3 ik lm - k lk kuk 3m ] C^ 11)} (p, k 3 - p, 
775-^3 [k 2k k 2 ik 3m - k 3k k 3l k 2m ] C^ n)l (p, k 2 -p, 
[hkhihm - kikkiik 2m ] C^ n)l (p, k 2 -p, 
3 [ki k kuk 3m - k 3k k 3 ik lm ] c|^ 11)} (p, k ± - p, 



(D.5) 



(2vr) 3 

d 3 p 

(2^) 



/ (Tl 3 [ fclfcfc i^2m - k 2k k 2 ik lm ] CH* (P> k x-p, 



-k 3 ) 
-k 3 ) 
-k 2 ) 
-k 2 ) 
-fci) 
-fci), 



2:? 



1 



-Bll®ll®ll — n {klmkln + k2 m k2n + k 3m k 3n/ i ,„ x 



3 ^mn yy) 



.(11), 



,(11) 



(11)/ 



,(11), 



kiik 2j c\- L >{k2) hkhtC^'iks) + fciifej- Q. li; (fci) fc 2fc %qr(fc 3 ) 



,(11) 



(11), 



+ fci*%^ Wfefc^q; 1 ^ 



+ 



+ 



+ 



-kuk 3j c[] 1] {k 3 ) j j^k lk k2iCf l 1 \p)k lm k2nC^ ) (p + k 2 ) 
^k 2i k 3j c\f ] \k 3 ) / T^3fc u .A;2/C , ^ 1) (p)/ci m A;2„C^ ) (p + fci; 
^fciifey C^Cfe) | 7& fc i* fc 3i C2 1} (p) k lm k 3n C£$(P + ks) 



+ ^iik 2j C^\k 1 ) 



y v 

,(n), 



(2vr) 3 
+ ^k 2i k 3j C^'(k 2 )l |^ 



k 2 fc fc 3 Z C'^ 1) (p) /c 2m A:3 n C^ (p + fc 3 



+ ^iifc 3 ,-c£ 1) (fci) 

d 3 p 



d 3 p 

(2^)3 



^lfcfc3Z Cklip) hmfon C^ (P + fc l^ 



k 2k k 3L C\^\p) k2mk 3n C^l\p + k 2 ) 



+ / j7j-^ 3 hik2j C^ n) (p) A;ife% Cm (P - fc i) fc2mfe3n c£V(p + fe 2 ) 



(D.6) 



Note that we have employed a shorthand notation, e.g., C 



{(112201113)} _ .-,{(1122)} ^,{(1113)} 



ijkl 



= a 



ijkl 



+CK 



ijkl ' 



where {•••} indicates summation over all possible correlator settings, e.g., C kl 
<C + C'S 1 ' + C%r + CgT • -d C<}P> is understood to mean IMI C«f"> . A, 



usual, summation over repeated indices € {i,j,k,l} is implied. 

D.2 Non-Gaussian initial conditions 

The one-loop terms arising from non-Gaussian initial conditions are 



i&i^ijk 



D 



1120122! 



3113 = -ikuk 2j k 3k Cl- k 



{(11201220113)} 



(fel,fe 2 ,fe 3 



(D.7) 



5ii0i2 = k u k 2 jk 2l k 3m [cg 1} (fei) C}£ 2)} (k 3 ) + C« 12)} (fci) C\l]\k 

' ,(n) 



+ kuhjhih,,,, C^\k 2 ) C}i 12)} (k 3 ) + Cf 2)} (k 2 ) d£\hs) 
cf , ik 1 )C^(k2) + C^(k l) C^ik2) 



+ kuk 3i k 2 ik 



li«-3i«-2/«'3m 



(D.8) 



D 1 / 7 1 7 f ^P ^,{(111101112)}/ , I I \ 

#111101112 = -jkiikijk2kk 3 i I j^-^ C>> kl (p, fei - p, k 2 , k 3 ) 

+ -k2ik 2 jk lk k 3 i / —-3 C^ niffilll2)} (P. fc 2 - P, fei, fes) 

, 1, , , , /" d3 P .-,{(111101112)}/ , 1 1 \ 

+ ^k 3 ik 3j k lk k 2 i / ^-^3 C^ v w ;/ (p, &3 - P, fei, fe2j 



(D.9) 
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£n®m = #110211(11 ® 211 -»■ 11 ® 111) 

/j3 
7^X3 [fcljfclm + &2^2m + fc 3 /fc 3 m] C ta (?) 

+ i / j^-Tzk u k2j C^ n) (p) k lk k 2 ik 3m CJ^'fa -p,p + k 2 , k 3 ) 

+ i / J^)3 klik 3J C ^ l) (P) k lkk 2 lk 3m C^m fa ~ P> fc 2,P + &3) 

+ i y 703^% Cg 1} (p) fcifcfcaifcam cg^Cfei, k 2 -p,p + fc 3 ) (D.10) 

+ xfattaj (^-(fci) / ——3 [k 3k k 3l k 2m ~ k 2 kk 2 lk3m] C k ]m (P> fc 3 - P, -fe 3 ) 
+ ^lifo? C , f j L1) (fc 2 ) / 7^3 [hkhikim - k lk kuk 3m ] CJ^'(p, k 3 - p, -k 3 ) 
+ ^kuk 3j Clf-'fa) I ——3 [k2kk 2 ik 3m - k 3k k 3 ik 2m ] CJ^'(p, k 2 - p, -k 2 ) 

+ ^kuk 3j C>j >fa) / ——3 [k 2k k 2 ik lm - ki k k u k 2m \ C^ m '(p, k 2 - p, -k 2 ) 

/i 3 
——3 [k lk kuk 3m - k 3k k 3l ki m ] CJ. lm '(p, fci - p, -fei) 

+ rcfeifoj C^ n) (fc 3 ) / ——3 [k lk knk 2m - k 2k k 2i k Xm \ C^\p, k x - p, -k x ) , 

{(122)} {(122)} 

again with the understanding that C^A. = i©nCJA. , and so forth. 

E Diagrams and the relation between the SPT and the LPT bispectra 

In this section we demonstrate how to construct diagrams in LPT. Consider the diagram in 
figure 4, corresponding to the contribution: 

* 1 r ni ll )ti \l 1 1 /o(H2)/ , , seq.(B.4) 

- -kuk 2 jCy '(fci)£:3fc%«2mC^ m ; (p,fc3 -p, -k 3 ) = 

- kuSV (k^PUk^k^fa) 

x k 3k S k 1) (p)P L (p)k 3l sl 1) (k 3 -p)P L (\k 3 -p\)k 2m S^(p,k 3 -p) , (E.l) 

where for simplicity we have omitted the integration over p. The RHS of equation (E.l) can be 
split into five parts: (i) ^(feijP^^i), (ii) k^S^^P^fa), (iii) hkS^^P^ip), 

(iv) %^ (1) (fc 3 -p)P^ /2 (|fc3 -p|), and (v) k 2m S$(p,k 3 -p)P 1 L / \p)P 1 L /2 (\k 3 -p\). 

Part (i) corresponds to an external vector k\ entering at the point i, from whence 
emerges an internal line of momentum k\ (represented by S\ (k\)), terminating at the 
linear power spectrum Pi{k\). Parts (ii) to (iv) bear similar interpretations. Part (v) has 

an external line k 2 entering at the point m, and two internal lines with momenta p and 

(2) 
k 3 — p emanating from it (indicated by Sm (p, k 3 — p)), which terminate respectively at the 

power spectra Pl(p) and Pl(|A; 3 — p\). Parts (i) and (ii) are joined together at their common 

linear power spectrum Pl(/ci), thereby linking the point i to the point j. Likewise, parts (iii) 
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U k 




Figure 4. An example diagram contributing to the one-loop bispectrum. The letters {?, j, /c,/,m} 
are indices of the correlators CV,-fc..., while ki + k 2 + k 3 = are the wavevectors defining the triangle 
of the bispectrum. The symbol © denotes an initial linear power spectrum. 

Table 1. The number of LPT contributions to each class of SPT diagrams for Gaussian initial 
conditions. 
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-6110211 
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B222 
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I-B321 


18 
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18 


6 


II-B321 


12 


6 


18 


6 


-B411 


12 


6 


15 


3 



is joined to one of the internal lines of part (v) at their common power spectrum Pl(jj), 
thus connecting the points m and k. The joining of parts (iv) and (v) at the power spectrum 
-Pl(|&3 — p) results in a connection between the points m and I. Finally, the points k and I are 
fused into one single point because of the common external vector k% in parts (iii) and (iv) , in 
the same manner that j and m are fused together because of the common external vector ki 
appearing in parts (ii) and (v). The resulting diagram contains a linear power spectrum in 
every internal line, and belongs to the class of 1-B321: three internal lines emerge from {j, m}, 
two from the point {A;, /}, and one from the point i (see also diagram (b) in figure 1). 

The construction of all other diagrams is then completely straightforward. Only one 
new rule needs to be introduced in the case of primordial non-Gaussianity: while the linear 
power spectrum has only one ingoing and one outgoing momentum, the initial bispectrum 
(trispectrum) has three (four) entries/exits. 

One advantage of using these quasi- Feynman rules is that it is now possible to associate 
all LPT contributions in equation (3.11) to their counterparts in SPT, since both approaches 
should lead to the same diagrams for the iV-point function. Note, however, that for every 
SPT diagram there are several LPT contributions. The number of LPT contributions to each 
diagram is listed in tables 1 and 2. 

Working with the diagrams in figure 1, we find the following correspondence between 
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Table 2. Same as table 1, but for additional contributions arising from non-Gaussian initial condi- 
tions. 
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the SPT kernels and the LPT correlators. For Gaussian initial conditions: 

B 222 {k u k 2 , fc 3 ) = -ifeCf 222 ^!, k 2 , k 3 ) + l -k ll2Z C^ ll22 \p, k ± - p, k 2 , k 3 ) 

+ -k 22 i 3 C (1122) (p, k 2 -p,k u fc 3 ) + -k 33 i 2 C (1122) (p, k 3 -p,k u k 2 ) 

+ \k 12 C^ l \p) fci 23 C( 112 )(fei - p, P + k 2 , k 3 ) 
+ ik 13 C^ n \p) fci 23 C (121) (fci - p, k 2 ,p + k 3 ) 
+ ik 23 cW ( P ) k 123 C^ (kx ,k 2 -p,p + k 3 ) 

- k l2 C^\p) k 13 C^(p - fei) k 23 C^(p + k 2 ) , (E.2) 



lSi23(fei,fe 2 ,fc3) = -ifci23lC (123) (fcl,fc 2 ,fc 3 



:(2)/ 



7 (l) 



W/ 



:( 2 ), 



k 33 i 2 S^(k h k 2 - p) S^(-p) SjW(fei) Sf\p, k 2 - p) 



?( 2 )/ 



:(1) 



(i) ( 



:( 2 )/ 



+ A:3312^ Z; (fcl,p)^ 1J (p-fc2)^(fcl)^ ; (p,fc2-p) 

+ fcasuSf >(-p) 5J 2) (-fei,p - fez) 4 1} (fei) ^/ (2) (-P,P - fea) 

+ fc3312<Sf } (P - *&) ^(-fei, -p) ^(fcl) 5f } (-p,p- fe 2 ) 
xPL(fcl)P L (|fc 2 -p|)P L (p) 

+ ^ 2213 C( 1113 )(p, fc 2 - p, k u fc 3 ) + fci 3 C (11) (fci) fe 23 C( 21 )(fc 2 
+ ^ 23 C (11) (P) A;i23C (112) (fei, fe 2 - P,P + k 3 ) 
+ \k 23 C^ ll \k 2 - p) fci23C (112) (fei,P, -P - fei) 
+ ^fciaC^^Cfei) [fc 223 - fc 332 ] C( 112 )(p, fc 2 - p, -fe 2 ) 
- ^ 13 C (11) (fei) k 23 C^ ll \p) k 23 C^ l \k 2 - p) , 



(E.3) 
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Il5l23(fcl,fc2,fc 3 ) = -ifcl23IlC (312) (fcl,fc 2 ,fc 3 ) 

-2k ll23 sf\-p, k x ) sf\p) s£\k 2 ) 5, (2) (-fei, -k 2 )P L (h) P L (k 2 ) Pl(p) 
-2fcii23^ (1) (p) 5J 2) (p, -fei) S^fe) 5f 2) (-fe ls -k 2 )P L (h) P L (k 2 ) P L ( P ) 

+k 13 C^ 1 \k 1 ) k 23 C^ n \k 2 ) + ifci 23 C( 112 )(fci, fc 2 , fc 3 ) fciiC( u >(p) 

+A ;23 C( 11 )(fc 2 ) fc 113 S< 2 } (-P, -fei) <S; (1) (-P) ^(-fci) ^l(p) ^(*i) 
+A; 23 C( 11 )(fc 2 ) fcusS^Cp) S, (2) (p, -fei) sW(-fci) Pl(p) P L (h) 
-ifciiC^^fcisC^Cfci)^^ 11 )^), (E.4) 

P 4 ii(fci, fc 2 , k 3 ) = -ik 123 C^ u \k h k 2 , k 3 ) + k 3312 C^ lll \p, k 3 - p, ki,k 2 ) 

-4*33i25? } (-fci,P) Sf\k 2 ,p) s£\kt) 5, (1) (fc 2 ) P L (p) P L (k!) P L {k 2 ) 
+k l3 C (ll \k l ) k 23 C^ Vi \k 2 ) + ki 3 C^ 13 \kx) feaC* 11 ^) 

+ ^i2 3 C (112) (fci,fc 2 ,fc 3 )A :33 C( 11 )(p) 

-2A; 13 C( 11 )(fc 1 ) fc 332 5f (p, -fc 2 ) 5, (1) (p) 5«(-fc 2 )P L (p) P L (fc 2 ) 
-2A; 23 C( 11 )(fe 2 ) fc33i5i 2) (p, -fei) ^(p) 5&>(-fci)P L (p) P L (ki) 
-\k 33 C^(p) fciaC^Cfex) fc 23 C (11) (fc 2 ) , (E.5) 

where we have employed the shorthand notation ki 23 C^ 222 ' = kuk 2 jk 3 kC\- k and so on, and 
omitted writing out the integration over p. Note that these expressions represent only one 
permutation of the full solution, e.g., for the full iP 32 i one must permute the wavevectors in 
l-Bi 2 3(fci ; k 2 , k 3 ) five times and then sum the six terms together. The kernels S^ a ' contain in 
general both longitudinal and transverse components. Contributions to B from transverse 
components vanish if and only if p = p\ 2 ... a in p ■ S^ a '(pi, . . . ,p a ). As it turns out, the only 
non- vanishing transverse contribution appears in the second term in equation (E.5). 
For non-Gaussian initial conditions we find: 

nBii2(ki,k 2 , fc 3 ) = -ifcira nC (112) (fci, fe, fe 3 ) + \k 33l2 C {llll \p, k 3 - p, fei, fc 2 ) , (E.6) 

lPl22(fel, k 2 , fe 3 ) = -ifcl 23 lC( 122 )(fc!, fc 2; fe 3 ) + kuC^ik!) k 23 C^(k 3 ) 

- k 3312 sf\p) sf >(fc 3 - p) S%\k x ) P L (h) 5f } (fei, fe 3 ) B (k 3 ,p, \k 3 - p\) 

+ ^i2C( n )(fei) k 332 C^ lll \p, k 3 - p, -fe 3 ) , (E.7) 

IlSi22(fel, &2, k 3 ) = ~ifel23 IlC (122) (fei, fc 2 , fe3) 

~2k 2213 S { i l) {k 2 -p)P L (\k 2 -p\)sf ) {p)S ( £ ) {k 1 )S ( i 2 \k 2 -p,p + k l ) 

2fc 33 i 2 5; i) (p - fc 2 ) P L (\k 2 - P\) 5 J (1) (-fci - P) 5i 1) (fei) Sl 2) (p - k 2 , -p) 
k 23 C^\k 2 -p) k 123 S i ^\k l )s\ 1 \p)S^(-p - fei)l B (h, P , |fci +p|) , (E.8) 
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iB 11 3{ki,k 2 ,k 3 ) = -ifci 23 iC (311) (fci,fe 2 ,fc 3 ) - h ll C^ ll \p)B {k l MM 

- k ll2Z S { p{p) Sf\p, -k x ) s£\k 2 ) 5 ; (1) (%) P L (p) B (h, fe, fe 3 ) 

- k ll23 sf\p, fci) Sf\p) S ( k 1) (k 2 ) ^ (1) (fe 3 ) Pl(p) B (h,k 2 , fe 3 ) , (E.9) 

nSnaCfci, k 2 , k 3 ) = -ife 123 n C( 113 )(fc 1; fc 2 , fe 3 ) + fci 3 C (11) (fci) fe 23 C( 12 )(fc 2 ) 

- fe 3 3i 2 5f ^fci, k 2 - p) Sf\-p) S^\k{) P L (ki) ^ (1) (fe 2 ) S (P, fe, |fc 2 -P|) 

- fe 3 3i2^ (1) (p - k 2 ) Sf\k u p) SJpfa) P L (h) ^ (1) (fc 2 ) B (p, k 2 , \k 2 - p\) 

- l -k l3 C^{k x ) k 332 C^ ul \p, k 2 - p, -k 2 ) . (E.10) 

Again, we emphasise that these expressions need to be summed over all their permutations. 
In particular, the summation over the three permutations of the term n-Bn 2 (fei, fe 2 ,fe 3 ) is 
missing in the current SPT literature [18, 32]. 

F Density kernels from LPT 

The symmetrised kernels for equation (3.18) are 

X[ s \k- Pl ) = k-S^(pi), (F.l) 

X { 2 s \k; Pl ,p 2 ) = k ■ S( 2 \ Pl ,p 2 ) + \k ■ S«(pi) k ■ SW( P2 ) , (F.2) 

X { 3 s \k; Pl ,p 2 ,p 3 ) = k ■ S^( Pl ,p 2 ,p 3 ) + \k ■ S«(pi) k ■ S^(p 2 ) k ■ SU( P3 ) 

b 

+ -{fe-5 (1) (pi)fc-5 (2) (p 2 ,p 3 ) + two perms.} , (F.3) 

Xi s \k;p 1 ,p 2 ,p 3 ,p 4 ) = k • S ,(4) (pi,p 2 ,p 3 ,p 4 ) 

+ ^k ■ S^ipi) k • S«(p 2 ) k ■ S^(p 3 ) k ■ S«(p 4 ) 

+ 3 I 2 fc ' S(2)( -P^P^ k ' ^ (2) (P3,P4) + two perms. I 
+ -< k- S ,(1) (pi)fc-S ,(3) (p 2 ,p3,p 4 ) + three perms. > 

+ I 1 \ k ' 5(1)(pi) k ' SW( -P^ k • 5(2) (P3,P4) + five perms. I . (F.4) 

At this point it is important to ask when purely transverse Lagrangian modes affect the nth 
order density contrast. The first term on the RHS of equation (F.3) is k ■ S^ 3 '(pi,p 2 ,p 3 ), 
where the S^ kernels contains a transverse term oc u> 3c (p\,p 2 ,p 3 ) (see equation (2.13)). 
At third order, the Dirac delta in equation (3.18) fixes k = p\ 23 , thus leading to pi 23 • 
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w 3c (P1jP2jP3) = because of the transverseness condition (A. 12), i.e., the third order 
transverse kernel u^J does not contribute to 5^'. However, when a similar k- S^'(pi 1 p2,Ps) 
term appears in the 5^' expression (F.4), the Dirac delta now implies k = P1234, so that 
P1234 • t*4c (Pi>P2,P3) = Pi • ^3c (Pi)P2,P3)) which is generally non-vanishing. Hence, the 
transverse term u^J does contribute to 5^ 4 ', and we can conclude by induction that all nth 
order transverse kernels must contribute to §( N>n > , but not S^ - n '. 
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